// BezierExtrem.cpp : Defines the entry point for the console application.
//

#include "stdafx.h"
#include <stdlib.h>
#include <time.h>
#include <float.h>

#define  MAX_DEG	21


void deCasteljau(double u, int n, double* func, double* re)
{
	double u1 = 1-u;
	int o=0;
	double nn[MAX_DEG];
	if (n>=MAX_DEG-1)
	{
		printf("Error!\n");
	}

	for (int k=0; k<=n; k++)
			nn[k] = func[o++];
	for (int k=1; k<=n; k++)
		for(int l=0; l<=n-k; l++)
			nn[l] = u1*nn[l] + u*nn[l+1];


	*re = nn[0];

}

void deCasteljau3(double u, double v, double w, int nx, int ny, int nz, double* func, double* re)
{
	double nn[4];
	double f1[4][4];
	double f2[4];

	double u1 = 1-w;
	int o=0;
	for (int i=0; i<=nx; i++)
	for (int j=0; j<=ny; j++)
	{
		for (int k=0; k<=nz; k++)
			nn[k] = func[o++];
		for (int k=1; k<=nz; k++)
			for(int l=0; l<=nz-k; l++)
				nn[l] = u1*nn[l] + w*nn[l+1];

		f1[i][j] = nn[0];
	}

	//GetBernstein(nn, ny, v);
	u1 = 1-v;
	for (int i=0; i<=nx; i++)
	{
		for (int k=0; k<=ny; k++)
			nn[k] = f1[i][k];
		for (int k=1; k<=ny; k++)
			for(int l=0; l<=ny-k; l++)
				nn[l] = u1*nn[l] + v*nn[l+1];

		f2[i] = nn[0];
	}

	u1 = 1-u;
	for (int k=0; k<=nx; k++)
		nn[k] = f2[k];
	for (int k=1; k<=nx; k++)
		for(int l=0; l<=nx-k; l++)
			nn[l] = u1*nn[l] + u*nn[l+1];

	*re = nn[0];

}

int _tmain(int argc, _TCHAR* argv[])
{
	int maxDeg=MAX_DEG-1;
	int maxExa=1000;
	int part = 100;
	int PartNum;
	double re, maxP, minP, maxPart, minPart;
	int maxPI, minPI, maxPartI, minPartI;

	srand((unsigned)time(NULL));

	for (int deg=2; deg<maxDeg; deg++)
	{
		double *p = new double[deg+1];
		int exa = 0;
		for (; exa<maxExa; exa++)
		{
			maxP = -DBL_MAX;
			minP = DBL_MAX;
			minPI = maxPI = -1;
			//for (int i=0; i<=deg; i++)
			//{
			//	p[i] = 2.0*((double) rand()/(double) RAND_MAX)-1.0;
			//}
			//p[0] = -0.63908810693685725;
			//p[1] = -0.63481551561021754;
			//p[2] = -0.68114261299478129;
				
			for (int i=0; i<=deg; i++)
			{
				p[i] = 2.0*((double) rand()/(double) RAND_MAX)-1.0;
				if (p[i]>maxP)
				{
					maxPI = i;
					maxP = p[i];
				}
				if(p[i]<minP)
				{
					minPI = i;
					minP = p[i];
				}
			}
			if (maxPI==-1||minPI == -1)
			{
				printf("Error in P\n");
				exit(0);
			}

			PartNum = deg*part;
			maxPart = -DBL_MAX;
			minPart = DBL_MAX;
			maxPartI = minPartI = -1;
			for (int parti=0; parti<=PartNum; parti++)
			{
				deCasteljau(double(parti)/PartNum, deg, p, &re);
				//printf("%d : %f\n", parti, re);
				if (re>maxPart)
				{
					maxPart = re;
					maxPartI = parti;
				}
				if(re<minPart)
				{
					minPart = re;
					minPartI = parti;
				}
			}
			if(maxPartI==-1||minPartI==-1)
			{
				printf("Error in parti\n");
				exit(0);
			}
			if (abs(maxPartI-maxPI*part)>part||abs(minPartI-minPI*part)>part)
			{
				printf("Error in deg = %d \n", deg);
				break;
			}
		}
		if (exa==maxExa)
		{
			printf("OK in deg = %d \n",deg);
		}
		delete []p;
		p = NULL;
	}

	return 0;
}

